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I. INTRODUCTION 

Superconduction correlations of electrons in weakly 
disordered metallic nano-grains^'^ are quite generally de- 
scribed by the Richardson Hamiltonian. This model was 
first introduced in the sixties in the context of nuclear 
physics'^''* and is Bethe Ansatz solvable. This means that 
the spectrum is known in terms of a set of parameters 
whose possible values are determined as the simultane- 
ous roots of a large system of (non-linear) equations that 
have to be solved numerically. 

Most studies of superconduction correlations in nano- 
grains therefore heavily rely on numerical solutions of the 
Bethe Ansatz equations for the Richardson Hamiltonian. 
What remained lacking was a simple intuitive picture of 
the low-energy physics of such grains, that would come 
from analytical rather than numerical studies. This was 
addressed by Yuzbashyan, Baytin and Altshuler^ who de- 
veloped a strong coupling expansion based on the Bethe 
Ansatz solvability of the Richardson Hamiltonian. Thus 
the authors were able to compute the spectrum in orders 
of the the inverse of the pairing interaction strength. 

In this paper we complement the work ofRef.iby ex- 
panding around the same limit but without relying on the 
Bethe Ansatz solvability of the Hamiltonian. As a result, 
perturbations that destroy the Bethe Ansatz solvability 
of the model can also be considered. Our strategy is 
to perform a so-called generalized Dyson boson-fermion 
mapping on the Hamiltonian. Generalized Dyson boson- 
fermion mappings jSiSiSiiS or Dyson mappings for short, 
are instances of boson expansion methods.^ A Dyson 



mapping is an invertible linear transformation that maps 
a many-fermion system onto a system in which indepen- 
dent bosons and fermions are present. The bosons in the 
mapped system fulfill the role played by certain fermion 
pairs in the original system. For the Richardson Hamilto- 
nian, zero-momentum spin singlet Cooper pairs are thus 
represented by bosons. 

The paper has three main aims: 

1. We want to extend the work of Yuzbashyan et al 
and build further on the intuitive picture for the 
low energy physics of superconducting nanograins 
by considering perturbations such as an external 
electric field that breaks the Bethe Ansatz solvabil- 
ity. We consider not only the spectrum, but also 
perform (zero-temperature) linear response calcu- 
lations for time-dependent perturbations. 

2. Superconductivity is often thought of as the Bose- 
Einstein condensation of Cooper pairs. However, 
Cooper pairs are not in themselves "fundamental" 
particles and can fairly easily be broken up into 
their constituent electrons by external perturba- 
tions. Viewing the superconducting state as a con- 
densate of non- or weakly interacting bosons leads 
one into the temptation to suppose that the break- 
up of a Cooper pair is an event that does not de- 
pend on the state of the rest of the system. We 
point out that this picture is wrong due to the ex- 
clusion principle. Because Cooper pairs consist of 
fermions, a high density of Cooper pairs reduces 
the number of unoccupied single particle orbitals 
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into which the electrons resulting from Cooper pair 
break-up can decay. As a result, the linear response 
of the system to a pair-breaking perturbation is not 
extensive in the number of Cooper pairs. 

3. We hope to convince readers of the utility and ele- 
gance of the Dyson Mapping as a tool for analysing 
the Richardson Model. To this end we uncover a 
selection rule for the matrix elements of single par- 
ticle operators between eigenstates of the strong- 
pairing limit using the mapping. This selection 
rule seems to have gone unnoticed until now, posi- 
bly contributing to the assertion in Ref. that a 
conventional perturbation analysis (not exploiting 
Bethe Ansatz solvability) is not feasible. Also, we 
show how the mapping simplifies the task of calcu- 
lating certain non-zero matrix elements of one-body 
operators between correlated many-body states. 

The paper is structured as follows. Section ^ takes 
care of preliminaries. The Richardson Hamiltonian and 
the Dyson mapping are introduced. The Hamiltonian 
obtained from the Dyson mapping when transforming 
(or mapping) the Richardson Hamiltonian is presented. 
In Section IIIII we motivate our claim that the Dyson 
mapping facilitates an analysis of the Richardson Hamil- 
tonian when the pairing term dominates the one-body 
term. In Section IIVI a perturbative diagonalization in 
orders of the kinetic energy is discussed. We calculate 
the first three corrections to the ground state energy, 
and compare the results to those previously obtained!^ 
We also calculate the ground state expectation value of 
an arbitrary additive one-body operator to first order in 
the kinetic energy. Brief mention is made of results ob- 
tained for excited states. Section El is devoted to the 
study of time-dependent phenomena. The main calcu- 
lation in this section yields the linear response of the 
expectation value of an arbitrary additive one-body op- 
erator to an external perturbation. An appendix contains 
a detailed exposition of Rayleigh-Schrodinger perturba- 
tion theory for operators whose unperturbed spectra are 
degenerate. The formalism is general enough to hold for 
non-hermitian operators and is presented in a way that 
makes it clear how to proceed in principle to arbitrarily 
high order corrections. 

Surprisingly, we could not locate a ready source in the 
existing literature that presents degenerate perturbation 
theory in sufficient detail and general enough for present 
purposes. We have therefore included a more complete 
discussion in the appendix than is strictly necessary for 
the derivations in the main body of the paper. 

II. HAMILTONIAN AND MAPPING 

With the single particle states, l^cr), = + 
1, . . . , j — 1, j; a — ±, we associate fermion creation and 
annihilation operators a^^. and a^o-, that obey the usual 
anti-commutation relations. The index a distinguishes 



between spin up and spin down, while refers to mo- 
mentum or some other quantum number, such that the 
state \^J-,cr) transforms to |— /i — cr) under time-reversal. 
The operators and a^cr act on fermion Fock space 
where we denote the vacuum by |0). 

The Richardson Hamiltonian is constructed in the 
fermion space as follows. First an operator 

(1) 

is defined, with — 2j + 1 the number of single-particle 
orbitals. This choice normalizes |0) to unity. With the 
aid of the operator and its conjugate S, the Richard- 
son Hamiltonian can be stated as 

3 

H= J2 (aj+a.+ + al,_a_._) - GS+S. (2) 

It is time-reversal invariant and contains a one-body ki- 
netic term and a two-body interaction or pairing term. 

The Hamiltonian describes electron-electron interac- 
tions in weakly disordered metallic nano-grainsi^ For 
simplicity, we have omitted the spin-exchange interac- 
tion, that should in principle be present as well; we 
are primarily interested in superconduction correlations. 
Since the spin-exchange term commutes with all other 
terms in the Hamiltonian, it is straight-forward to re- 
introduce it at a later stage if necessary. Its ommision 
here does not represent a serious loss of generality. 

The following mathematical features of the model are 
important. 



A. The Seniority Limit 

If the kinetic term is switched off, the Richardson 
model reduces to Racah's seniority model, so that we 
call this the seniority limit. The spectrum is explicitly 
known in the seniority limit ^^ii^ For a system with 27V 
particles, (with N an integer), the eigen-energies coalesce 
into -I- 1 levels 

Ei°^ = -^iN-^m + ^-N-s), (3) 

with s = 0, 1, . . . , A. Here the superscript is intended 
to indicate that this is the spectrum in the absence of 
the kinetic term. The ground-state, with s = and en- 
ergy £;^'" = -GN{n + 1- N)/nis non-degenerate. All 
other levels are highly degenerate. The energy difference 
(known as the superconducting gap) between the ground 
state and the (degenerate) s = 1 level of eigenstates is 
G. Up to normalization, the ground-state eigenvector is 

{S+f\0). 
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B. Bethe-ansatz Solvability 

It is possible to diagonalize the Richardson Hamilto- 
nian beyond the seniority Hmit by exploiting the fact that 
it is Bethe-ansatz solvable. See for instance Appendix 
B of Ref. H We focus on diagonalization in the sub- 
space of time-reversal invariant states. These are states 
for which no single particle state (^fi) is occupied with- 
out the time-reversed state (— /i — cr) also being occupied. 
It can easily be shown that the Richardson Hamiltonian 
leaves this space invariant, and also that, if wc know 
how to diagonalize the Richardson Hamiltonian in the 
time-reversal invariant subspace, the diagonalization in 
the whole fermion space can be achieved without further 
complication, thanks to the so-called blocking effect. 

As a result of its Bethe-anzatz solvability, the follow- 
ing holds for the eigenstates and eigen-energies of the 
Richardson Hamiltonian in the time-reversal invariant 
subspace: for a system of 2N fermions, let {Ek}^^^ 
be a set of complex numbers that satisfy the so-called 
Richardson equations 

G ^.Sfc-2e, +,4^^ Ei-Ek' ^ ' 

for = 1, 2, . . . , A^. It can be shown that there are (^) 
such sets one for each dimension of the time- 

reversal invariant subspace. To each solution corresponds 
an eigenvector of the Richardson Hamiltonian, given by 

up to normalization. The energy associated with this 
eigenstate is 

N 

E = Y.Ek. (6) 

fe=i 

Actually finding all the sets {E^}^^^ that solve the 
Richardson equations amounts to finding the simultane- 
ous roots of A^ polynomials of A^ variables and degree 
Q, + N . If either the number of single particle orbitals 
(fi) or the number of fermions (2A^) becomes too large, 
or if analytical rather than numerical results are required, 
the Richardson equations become an obstacle. 

C. Strong Coupling Expansion 

As a consequence, few explicit results have been de- 
rived for this model, and most studies involve numer- 
ical results based on the exact solution. To redress 
this state of affairs, Yusbashyan, Baitin and Altshuler 
(YBA)recently developed a perturbation expansion in in- 
verse powers of the coupling constant G. With this ex- 
pansion analytical results for the spectrum can be ob- 
tained away from the thermodynamic limit, A^, 1/ — > oo, 



where the discreteness of one particle levels becomes im- 
portant. In principle, the strong pairing limit, consid- 
ered by these authors, is realized in a small grain of ideal 
regular shape so that single electron levels are highly de- 
generate. In the case where the grain is so small that 
the energy distance between degenerate many-body lev- 
els is much larger than G, only the partially filled Fermi 
level is relevant. The kinetic term is simply a constant 
proportional to the total number of particles, and can 
be set equal to zero. Deviations from a perfectly regular 
shape slightly lifts the degeneracy of the Fermi level and 
is incorporated by treating the kinetic term as a pertur- 
bation. 

Of course, in an irregularly shaped grain, the kinetic 
term cannot be viewed as a small perturbation to the 
pairing term. However, even in this regime a pertur- 
bative analysis of the strong coupling limit provides an 
intuitive insight into the nature of superconducting cor- 
relations. Whereas the pairing term cannot be treated 
perturbatively, because the radius of convergence for the 
expansion shrinks to zero in the thermodynamic limit, 
the 1/G expansion has a finite convergence radius. Evi- 
dence suggests that for realistic 1/G the system behaves 
qualitatively similar to the small 1/G limit. 

In the next sections, we will do such a 1/G perturba- 
tive analysis of the model. Our expansion method com- 
plements that of Ref. 0, in that it does not rely on the 
Bethe-Ansatz diagonalization of the Hamiltonian. Thus, 
for instance, our method remains applicable when the 
kinetic term is replaced by an arbitrary one-body per- 
turbation, even though the resulting Hamiltonian can no 
longer be diagonahzed using the Bethe-Ansatz. 



D. Mapping the Building Blocks of the Model 

Our startegy involves mapping the Richardson Model 
onto an equivalent boson- fermion model, by means of the 
so-called Dyson mapping. The mapping uncovers some 
hidden features of the model and simplifies the task of 
calculating matrix elements of observables between cor- 
related many-body states. (See the Refs. HollTsI for a 
more detailed discussion of the properties of the Dyson 
mapping.) 

The operators S and 5+ are pseudo-spin ladder opera- 
tors for an SU{2) irreducible representation of dimension 
il, or in other words of total pseudo-spin As 
such, they have a lot in common with the boson creation 
and annihilation operators and B, especially if Q is 
large. The Dyson mapping is a linear invertible operator 
designed to replace the creation (annihilation) of Cooper 
pairs with the opperator 5*+ (5), with the creation (anni- 
hilation) of bosons. Without showing how the mapping 
is derived we simply give the (linear extensions of) the 
Dyson images of the components needed to construct a 
system described by the Richardson Hamiltonian. (See 
Ref. for a complete discussion on how the Dyson map- 
ping is defined as well as a derivation of the images of 
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the operators that constitute the building blocks of the 
I 



In these expressions, aj,^ and its conjugate are fermion 
operators, with the usual anti-commutation relations. 
We will refer to the the particles created and annihilated 
by them as ideal fermions, to distinguish them from the 
fermions of the unmapped system. and B are bo- 
son operators with the usual commutation relations. Im- 
portantly, the boson operators commute with the ideal 
fermion operators. Furthermore, counts the number 
of ideal fermions, Mb counts the number of bosons and 
is defined as the operator with the real fermion 
operators a+g. replaced by ideal fermion operators aj^^.. 
The fermion vacuum |0) maps onto the boson- fermion 
vacuum |0), while the functional (0|, maps onto a func- 
tional (0|. 

Note two subtleties of the mapping. Firstly it is not 
unitary, in the sense that (X)^ ^ (We de- 

note hermitian conjugation in the fermion space with a 
-|- sign while in the boson fermion space we denote it 
with a f, to emphasize that the mapping of the conju- 
gate is not the same as the conjugate of the mapping.) 
This implies that a hermitian fermion operator such as H 
will in general be mapped onto a non-hermitian operator 
{H)jy. The equivalence between the Richardson Hamil- 
tonian and the boson-fermion system is therefore a very 
formal mathematical one. It is not the physical equiva- 
lence between two systems that can both be realized in a 
(gedanken) laboratory. Secondly, the reader will have no- 
ticed the word "linear extension" used in connection with 
the Dyson images of fermion operators. The need for such 
terminology can be appreciated as follows: Recall first 
that the fermion Hilbert space on which the Richardson 
Hamiltonian is defined, was built up from a 2f7 = 4j + 2 
dimensional one-particle Hilbert space. The operators 
used in the construction of the Richardson Hamiltonian 
therefore act on a finite dimensional vector space. Yet the 
Dyson image of the fermion operator 5" for instance is the 
boson-operator B. While S acts in a finite vector space, 



Richardson Hamiltonian.) 



(7) 
(8) 



(9) 



(10) 

I 

B acts in an infinite-dimensional vector space: an infinite 
number of bosons can occupy the same state. What is 
happening here is that the operators in (|10|l leave a sub- 
space of the full infinite dimensional boson-fermion space 
invariant. Within this subspace the operators of H1(J|) 
are equivalent to the building blocks of the Richardson 
Hamiltonian. One might naively suppose that the physi- 
cal subspace of boson fermion space is that in which twice 
the number of bosons plus the number of ideal fermions 
equal the number of real fermions. However, counting di- 
mensions one sees that the dimension of this subspace of 
boson fermion space is still too large. The physical sub- 
space is still smaller and its identification more subtle, 
but fully solvedjlSi While this brings about the complica- 
tion that one should take care to choose initial and final 
states from this subspace, working in a larger than neces- 
sary Hilbert-space has advantages as well, as will become 
clear in the next section. 

At this point it is still not obvious what one gains by 
performing a Dyson mapping on the Richardson Hamil- 
tonian. This is addressed more fully in the next section. 
Note here, however, that the mapped operators have the 
following non-trivial but favorable properties: 

1. A state containing 2N real fermions maps onto a 
state for which twice the number of bosons plus the 
number of ideal fermions is equal to 2N. 

2. The operator {a^cr)^, only contains terms that anni- 
hilate an ideal fermion without changing the num- 
ber of bosons or that annihilate a boson and create 
an ideal fermion. In fact, with the annihilation of a 
real fermion is associated three coherent processes 
in the mapped system. In the first an ideal fermion 
is annihilated, in the second a boson is annihilated 
while an ideal fermion is created, and in the third 
an ideal fermion and a boson is annihilated while 
an ideal fermion pair is created. No terms that 
change the boson number by more than one occur. 



B, 



S^Ba^a 



n-Afp-AfB 
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Similarly, the operator (a^cr)^, '^^^ create an ideal 
fermion without changing the number of bosons or 
it can annihilate an ideal fermion and create a bo- 
son, but contains no terms that change the boson 
number by more than one. 

3. Also, a real fermion many-body state in which 
all fermions occur in time-reversed pairs a^_|_a^^„ 
maps onto a boson-fermion state in which 
ideal fermions only occur in time-reversed pairs 

4. A real fermion state for which the number of real 
spin up fermions minus the number of real spin 
down fermions is M, maps onto a boson-fermion 
state for which the number of ideal spin up fermions 
minus the number of ideal spin down fermions is 
also M. 

We stress that the extensions of the images of the oper- 
ators that appear above are exact and complete as they 
stand. The expressions are not the first few terms in se- 
ries of which 'higher order' terms have been neglected, 
nor was any other approximation made. 



Ad- 
Ado 

Ad+ 



In these expressions, the operator A is defined, similar 
to , as the operator A, with the real fermion operators 
a^g. and a^o- replaced by ideal fermion operators and 

a^ff . The operator A* = J2fj.,ya A-i,-^aj^^a,^a- is the op- 
erator into which A transforms under time-reversal. The 
operator is given by 

1 ^ 

and S,A is its hermitian conjugate in the mapped boson- 
fermion space Ti-BF- With triA we mean the trace 
^^j^ of the one-body matrix A^^. The separa- 

tion 01 Ad into three terms has been done in a way such 
that Ad- annihilates two ideal fermions and creates a 
boson while Ado leaves both the number of bosons and 
ideal fermions unchanged and Ad+ annihilates a boson 
and creates two ideal fermions. 



These building blocks can be combined to find the lin- 
ear extension Hd of the image of the Richardson Hamil- 
tonian. For the pairing term we simply get 

Pd - (S+S)^ = -^J^Bin + l-AfB- Nf). (11) 

Instead of immediately focusing on the image of the ki- 
netic term in the Richardson Hamiltonian, we rather con- 
sider an arbitrary additive one-body operator 

j 

A= A^, {a++a,+ + a+_a,-) , (12) 

that leaves the fermion space invariant^i and does not 
distinguish between spin up and spin down states. These 
two restrictions are not fundamental; they are imposed 
purely for the sake of simplicity. The linear extension of 
the image of A under the Dyson mapping can be conve- 
niently separated into three parts Ad = Ad- + Ado + 
Ad+, for which we have the following explicit expres- 
sions. 



(13) 
(14) 



(15) 



From H13|) - (|15|l it is clear that an additive one-body 
operator A = J^fiua- A^^a^^ai,fj that is odd under time- 
reversal, i.e. Afj_i, — — is simply mapped onto an 
operator Ad^ A^ Y.^,ua A tj.v0i^fia^f<y (Remember that 
hermitian operators that are odd under time-reversal are 
traceless.) This has the consequence that Ad commutes 
with Pd = —G{S~^S)jj. That additive one-body opera- 
tors that are odd under time-reversal commute with the 
pairing interaction, can also be checked without difficulty 
in the original system. 

To summarize the most important points of this sec- 
tion then, there is a similarity transformation Td , known 
as the Dyson mapping, that maps many-fermion states 
onto states in which Cooper pairs are replaced by bosons. 
Any fermion operator X transforms into an operator 
TdXTj^\ defined on a subspace of the full boson-fermion 
Hilbcrt space and is said to be the image of X. We have 
presented linear extensions of the images of the building 
blocks of the Richardson Hamiltonian to the whole boson 
fermion space. We denote the extended image of X un- 



A 



Mb 



{A + A, -2 tn A}~S^A 



B- 

{n - Mf - - jVf) 



{n - 2j\fB-AfF)i^ 
{n-AfF)i^-J^F + 2) 

{n-MB-MF)n 



-S^B{A + A, -2tri A) 



{n - jVf - - aTfY 



{n-jVB-MF)^^ 
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der the Dyson mapping as Xjj. Somewhat carelessly we 
often call Xo the image of X under the Dyson mapping. 
Bra's and kets in the original fermion Hilbert space are 
written with angular brackets, while bra's and kets in the 
new boson-fermion Hilbert space are written with round 
brackets. 

From the images of the operators in ir7|l- H10() the the im- 
age of the Richardson Hamiltonian under the Dyson map- 
ping can be constructed (|13|I - H15|) . In the next sec- 
tion we explain why it is useful to transform the Richard- 
son Hamiltonian using the Dyson mapping if the pairing 
term dominates the one-body term. 



III. THE UTILITY OF THE DYSON MAPPING 

Let us diagonalize the image Pd = — G (5*+ 5)^ of the 
pairing term, supposing that the original system con- 
tains 2N fermions. This implies that twice the num- 
ber of bosons plus the number of ideal fermions equal 
2A'' in the mapped system. Since the mapped pairing 
term is clearly diagonal in boson-fermion space, we 
immediately find by simple inspection that any state in 
boson-fermion space that contains N — s bosons and 2s 
ideal fermions is an eigenstate of Pq, with eigenvalue 
Ei°^ = -G{N - s){n + 1 - N - s)/n. This is ex- 
actly the same spectrum as that of the seniority model 
P = —GS^S on the fermion side, as is of course to be 
expected from a faithful mapping. However, if one counts 
the degeneracy of many-body levels, one finds that the 
mapped pairing term Pjj has a higher degeneracy than 
the original pairing term P. For instance the degeneracy 
of the first level of excitations (s = 1) for Pjj is — 51 
while for P it is 2U^ — 51 — 1. This discrepancy is due 
to the fact that Pd is defined in a larger Hilbert space 
than P. One of the eigenstates in the s = 1 level of the 
spectrum of Pd is unphysical. The identification of this 
unphysical state and of all unphysical states in the other 
levels of the spectrum of Pd has been carried out else- 
where p. 99]. The subspace of boson-fermion space 
on which the Dyson mapping establishes an equivalence 
with the Richardson Hamiltonian is completely known. 
For the purposes of the present text a full knowledge of 
this physical subspace is however not necessary. 



A. An Important Selection Rule 

We are now ready to state one of the main new tech- 
nical results of our investigation which subsequently also 
turns out to be of great utility. Every eigenstate (and 
therefore every physical eigenstate) of Pd with energy 

contains N — s bosons and 2s ideal fermions. As a 
remarkable consequence, the boson-fermion counterpart 
Ad of a, one-body operator A only has non-zero matrix 
elements connecting states in a level E'i"-' with states in 

P(0) 



Ad does not change the number of ideal fermions by more 
than two. Because of the linearity and invertibility of the 
Dyson mapping, this result is equally true in the original 
system: 



Ei°^ X E^f 







(17) 



unless s' G s — l;s;s + l, where \Es) is any eigenstate 
from level s of the pairing operator P = —GS~^S and 

similarly i^'^?''^ comes from level s'. X is any additive 
one-body fermion operator. Since this useful selection 
rule is not as immediately apparent in the original sys- 
tem, it seems to have gone unnoticed in the literature. 
With its aid the number of matrix-elements in a pertur- 
bation expansion for which non-trivial calculations are 
required is reduced by orders of magnitude. 



B. The Advantage of a Larger Hilbert Space 

When calculating matrix elements of operators, one 
often needs a spectral representation of the identity op- 
erator, or of a projection operator. (If {|fc)}fii is an 

orthogonal basis for a vector space then / = X^fcLi 1^) (^1 
is a spectral representation of the identity operator.) 

Now consider the somewhat artificial but helpful ex- 
ample of a three-dimensional vector space. Imagine that 
we are working in a 2-dimensional sub-space spanned 
by the two orthogonal vectors |1) — sin 9 sin(t)\x) + 
sin 6* cos(/) |y)-f cos6' \z) and |2) — cos(/) \x)~sin(j) \y). For 
vectors \a) and in this subspace, the inner product 
{a 1/3) can be written as (a|/|/3). The identity operator 
/ has the decomposition 



I = |1)(1| + |2)(2| 

= (sin^ 9 sin^ ^ -I- cos^ (/>) \x) {x\ 

— cos^ 9 sin (j) cos (f>\x) {y\ + seven similar terri^iS) 

Instead of this complicated expression, we are free to use 
the identity operator defined for the whole three dimen- 
sional space /extended = \x) {x\ + \y) {y\ + \z) {z\, since 
(a| / 1/3) = {a\ /extended |/3) for all |a) and |/3) in the two- 
dimensional subspace. 

This same advantage is to be had by working with 
the mapped instead of the original system: In any 
perturbation expansion where the pairing term is per- 
turbed by a one-body operator A, typically quantities like 
{ip's'l AIlsj^AIls2 ■ ■ - Ils^A \tps) must be calculated. Here 
Hs is the projection operator that singles out the com- 
ponent of the state it acts on in the space spanned by 
vectors that lie in level s of the spectrum of the unper- 
turbed Hamiltonian. The vectors \ips) and IV''^') denote 
eigenstates from the levels s and s' respectively. 

For simplicity, let us focus on the modestly non-trivial 
quantity {GS\ AIliA\GS) , which is encountered in the 
calculation of the second order correction to the ground 

N 



the levels i?^, with s' — s ~ 1, s, s + 1. This is because state in the strong pairing limit. Here \GS) cx {S+) |0) 
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is the normalized ground state in the seniority hmit. 
To make any head-way, an expression must be found 
for Hi. To do this, we need a complete basis for the 
s = 1 eigenspace. Finding such a basis is already a non- 
trivial task, and becomes progressively harder for larger 
values of s. It can be shown that any complete basis 
for the s = 1 eigenspace can be written in the form 



{\k) = {s+f-'j:i.__ 



k=l 



where toi 



is the degeneracy of the level s = 1 and Ckf^i, are complex 
numbers. (We only consider states with equal numbers 
of spin up and spin down fermions because the ground 
state contains as many fermions with spin up as spin 
down, and we have conveniently assumed A not to intro- 
duce spin flipping.) The projection operator Hi can then 
be written as IIi = Sfc=i l'^) (^1- This leaves us with the 
unpleasant task of calculating matrix elements such as 



{GS\A\k) = 



- |0) 



(19) 



In order to do this calculation, the operators must be 
normal-ordered. As far as we can see, the operators A, 
S, and "^^^ Ckfj.^O'^+CLt- do not belong to any useful 
algebra that would simplify the normal ordering. This 
means that, in the myriad commutation operations in- 
volved in normal ordering, new linear combinations of 
fermion operators that have not been encountered before 
are frequently generated. The result is that the num- 
ber of non-trivial operations involved in calculating l|19l) 
grows explosively as N increases. 

Let us now consider the calculation of 
{GS\AIliA\GS), using the Dyson mapping. Firstly, by 
replacing 5*+ with (5*+)^ and S with Sd we find that 



Td \GS) = C 



N 



\0)^C\-)., 



and 



{GS\T^'^C{0\ 



B 



N 



(20) 



(21) 



where C and C are complex numbers. The notation is 
intended to indicate that | — ) contains no ideal fermions. 
Because it holds that 1 = {GS\GS) = {GS\ T^^Td \GS), 
we know that CC — 1. Thus we find 



{GS\AniA\GS) 



lAn^A, 



(22) 



In this expression, the image of IliA\GS) is CAd+ |— )• 
The only part of {GS\ ATj^^ that has a non-zero overlap 
with this is C {GS\ Ad-- To compute (-| Ad-Ad+ |-), 
it is convenient to insert the identity operator in the 
s = 1 eigenspace (which is just a linear extension of 



J 



TqIIiTq^) between Ad- and Ad+- If we set |^i^*) 



(St) 



N- 



|0) then 



i-\AD-AD+\-) = 



1 

E 



{-\Ad- \^J.iy*) {^J-iy*\ Ad+ |-) . 

(23) 

Here the notation is intended to indicate that the 

state contains two ideal fermions (and therefore — 1 
bosons), and further that one of the ideal fermions occu- 
pies the state labeled by (/i, -I-) for which is shorthand, 
while the other occupies the state {—i', —) for which i/^ is 
shorthand. The matrix-elements in the above expression 
present no great challenge to compute. This is due to 
the fact that the bosons and ideal fermions have simpler 
commutation rules than the building blocks of states and 
operators on the fermion side. Particularly, ideal fermion 
operators commute with boson operators, as opposed to 
real fermion operators that have non-trivial commutators 
with the 5*+ and 5* operators. It is true that things be- 
come significantly more complicated in the mapped sys- 
tem for levels s higher up in the spectrum. However, 
there is no ceiling or clear boundary beyond which cal- 
culations suddenly become impossibly complicated. For 
simplicity we will only calculate quantities that involve 
matrix-elements between states from the levels s = and 
,s = 1 in this paper, but calculations involving states from 
the level s = 2 are no less feasible. Only, the calculations 
become tedious and the results too intricate for casual 
interpretation. 

We now write down all the matrix-elements that we 
are going to use in the rest of this paper. 



HAdoH = 2N 



tTiA 



(24) 
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{^\AD-\^iv.) = 

{iiv^\Ad+\~-) = 



N 

n 

^^-^ >A +A 



tri A 



TV - 1 

'2— :- tri A Sa'nSu'v + { 



N- 1 



n-2 



N ~ 1 

n-2 



n-2N 
'n{n - 2) 



(25) 
(26) 



{A-^f_l + A — ^ — jy) 



(27) 



Calculating the corresponding many-body matrix- 
elements in the unmapped system seems intractable to 
us. 



IV. TIME-INDEPENDENT PERTURBATION 
THEORY 



In this section we discuss the implementation of 
Rayleigh-Schrodinger (timc-indcpendcnt) perturbation 
theory for the linear extension Hd of the image of the 
Richardson Hamiltonian. The image Kd of the kinetic 
energy operator is treated as a perturbation to the im- 
age Pd of the pairing operator. The original kinetic en- 
ergy operator K defined on fermion space is an addi- 
tive one-body operator with one-body matrix elements 



Sfj.viii- It is assumed to be time- reversal invari- 



ant, I.e. e. 



We focus on the first few corrections 



to the ground state energy and eigenvector. Brief men- 
tion is made of the procedure involved in the expansion 
for the degenerate excited states. 

In a recent paper by Yuzbashyan and co-workers )^ a 
procedure was presented for determining the eigenval- 
ues perturbatively. The strategy there was to develop 
an expansion for X^feLi -^fe where Ek, k — 1, . . . , N sat- 
isfy the Richardson equations I0J. We do not claim that 
the method described in this section is easier to imple- 
ment than the method of Refj^ if one only wants to cal- 
culate corrections to the energies. What we do claim is 
that Rayleigh-Schrodinger perturbation theory is simpler 
to implement than the authors of Ref'F may have real- 
ized. This we attribute to the structure resulting from 
the Dyson mapping., as we proceed to show. With our 
procedure, it is feasible to calculate the first few cor- 
rections to eigenstates at the bottom of the spectrum. 
Thus, our procedure allows one to compute, for instance, 
the first few corrections to the ground state expectation 
values of observables. 



A. Ground State Energy 

As a demonstration that the formulas of the previous 
section are valid, we will start by calculating the first 
three corrections to the ground-state energy, and com- 
pare the results with those obtained by Yuzbashyan and 
co-workers. The general formulas for these corrections 
(from the appendix) are 



^(1) 

^GS 

^(2) 
^GS 

^(3) 
^GS 



HKdH, 

HKdQgsKdH, 



HKdQgs ( Kd - Ag^ ) QgsKd 



(28) 
(29) 
.(30) 



Here Qgs is a linear operator defined through its im- 
age on eigenstates of the mapped pairing operator Pd as 
follows. 



Qgs IV", s) 







if s = 
1-0, s) if s > 



(31) 



for any state {ip, s) with N — s bosons and 2s ideal 
fermions. Using the decomposition K]j = Kd- -\- Kdo + 
Ku+ together with this definition of Qgs-, we can replace 
(EHl-lEnil by 



^GS = i-\KDn\-), 



Ag) = {-\Ku- 



Kd+\-), 



^(3) 
^GS 



{-\Kd- [K 



^GS 



(32) 
(33) 

Kd+H- (34) 



(Remember that e["^ - E^°'^ = 1. 
from (j21|) that 



Immediately, we have 



2Ne, 



(35) 



where, in general, e'' — ^Ji/^- order to calcu- 

(2) (3) 

late Apg and A^,^, we insert appropriate representations 
of the identity operator to reduce the expressions in 133|) 
and (|34|l to sums over the matrix-elements in (|24|) - (|27|l . 
For the second order correction we find 



(36) 
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where the matrix elements are (— | if^^ — To calculate the third order correction, we need the 

2-\/iV/fie^ and |— ) = (e^ — e). The additional matrix elements (/i/x*|X£)o — ^^gs\^'^*) — 

fact that, in the sum over states, only states whose 2 "^^^ {(e^ - e) (5^^. - -ff- } . Thus we find 

ideal fermions form time-reversal invariant pairs appear, 
can be traced back to the time-reversal invariance of the 
original kinetic energy operator K. Hence we find 



= 7^ E {-\Ku-W.)(.mAKu^-'^%\vv:)(yv^KuA-) 



G2 

8A^(f7-iV)(f]-2iV) 
G217(17-l)(f]-2) 

8N{n~N){n~2N) ——^ 
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After converting the notation used in this paper to the 
notation of Ref. H namely N M , Q N , 7^ 
'2,~^Sk/N and G — *■ iVAd, we find our first and second 
order corrections agreeing with those in Ref. ^ The 
seemingly different third order corrections are the result 
of a small arithmetic or typing error in Ref. |^ It is 
not hard to check explicitly that the method of Ref. IH 
when implemented correctly, does indeed yield the same 
answer as we get in 



B. Grounstate Eigenvector and Expectation Values 

Let us now calculate the first order correction to the 
ground-state eigenvector. (To us it is not clear whether 
this can be calculated using the method of Ref. bmce 
the perturbation is not hermitian, a distinction must be 
made between the left and right perturbed eigenvectors. 
(Such a distinction was not necessary for the unperturbed 
ground-state as the unperturbed operator Pd happened 
to be hermitian.) The right ground-state eigenvector 
of the perturbed operator Ko + Pd is, to first order. 



\GSn) = I-) + 



GS'^j^^ I whereas the left ground state 



GS 



(1) 



also to first order. 



eigenvector is \GSl) = |— ) + 
Here \GSl) is defined to be the ground-state eigenvector 
of {Kd + PdY- 

According to the appendix, the first order correction 
for the right-eigenvector is 



GS 



(1) 



QgsKdH 
1 



G 



Kd+\-) 



= 2i 



G 



■M, (39) 



where the second line follows from the selection rule for 
the matrix-elements of the images of additive one-body 
operators explained in Section IlIII Similarly we have for 
the left eigenvector 



G5i'M = QgsKI, 



3 



(40) 



At this point we emphasise that the eigenvectors found 
above are physical. If one were for instance to find the 
ground state perturbatively to first order for the origi- 
nal fermion Hamiltonian and then apply the Dyson map- 
ping, one would obtain the above right-eigenvector. In 
general, implementation of exactly mapped operators ob- 
tained from the Dyson mapping do not per se introduce 
spurious components. It is only when the Hamiltonian 
is diagonalized in the complete boson (or boson-fermion) 
Fock space, rather than in the physical subspace, that 
spurious states may be obtained. While the physical sub- 
space may be obtained ab initio via the mapping of a 
basis in the original space, this would nullify the relative 
simplicity obtained when using a standard basis for the 
complete boson (or boson-fcrmion) space. Nevertheless, 
if the Hamiltonian is mapped exactly, as we do here, (i.e. 
without any truncation of degrees of freedom), the re- 
sulting eigenstates neatly separate into completely phys- 
ical states and states with spurious components, and a 
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number of ways exist to distinguish these (see e.g. Refs. 

For the Richardson Hamiltonian this is done in 
in a novel way - see Ref. fisl p. 99. The neat division be- 
tween completely physical states and states with spurious 
components still holds when the mapped Hamiltonian is 
diagonalized perturbatively in the kinetic term, as is ex- 
plained in detail in the last paragraph of the appendix. 

To calculate higher order corrections to these eigenvec- 
tors, we need to know matrix elements between higher ex- 
cited states. Specifically, for the second order corrections, 
matrix elements of the kinetic term between states in the 



s = 1 and 5 = 2 levels are required. We will however not 
proceed beyond the corrections already calculated. 

Let us calculate the ground state expectation value 
for a one-body additive operator A to first order 
in the perturbation. It holds that {GS\A\GS) 
C {GSl\ Ajj \GSr), where C is a normalization constant 
such that 1 = C{GSL\GS]i). Since it already holds that 
{GSl\GSr) = 1 + O (e^/G^), and we are interested in 
corrections of order e^/G, we can take C = 1. Thus we 
have 



{GS\A\GS) = HAdo\-) + {gs[ 



A 



D+ 



An- 



GS T 



G{VL-l) 



(41) 



where A = tri A/fl. 

We briefly discuss this novel result: Superconductiv- 
ity is often presented as the Bose-Einstein condensation 
of Cooper pairs. The ground state contains all electrons 
bound into Cooper pairs with all Cooper pairs in the 
same quantum mechanical state. The first term in H41|) 
confirms this picture: each Cooper pair contributes 2A 
to the expectation value. One could then easily suppose 
that the kinetic perturbation breaks up Cooper pairs, and 
that each Cooper pair can be broken up independently 
from all the others, so that the first order correction 
should also be extensive, i.e. proportional to the number 
of Cooper pairs in the ground state. However the first 
order correction is actually proportional to N{fl — N). 
The fundamental reason for this A^-dependence is block- 
ing through the exclusion principle. There are only 2Q 
states for the fermions making up the Cooper pairs in 
the original system to occupy. Therefore, when there 
are 2N electrons, all one-particle states are full, and no 
rearrangement of electrons is possible. In this case the 
perturbation has no effect. In the mapped system this 
"taking account of the exclusion principle" manifests it- 
self in the complicated images of the kinetic operator that 
introduce interactions of ideal fermions with bosons and 
with other ideal fermions. The deviation from bosonic 
behavior where the breakup of Cooper pairs is concerned 
will be encountered again in the remainder of this text 
where we will further comment on it. 

As a very simple application of 1)41(1 , we find the prob- 
ability Pp for the single particle orbital p to be occu- 
pied, by calculating the ground-state expectation value 
of (a^^ap+ -t- ap_ap_) /2. According to H41|) we have 

N 2N(fl~ N) , , 
As the kinetic term is switched on, single-particle orbitals 



with below-average kinetic energies are populated with a 
higher probability, in order to get a smaller kinetic energy 
contribution. In other words, the distribution in (|42|) is 
the result of a competition between the pairing term, 
that strives for a distribution where all single particle or- 
bitals are equally likely to be occupied, and the kinetic 
term, that would like to populate the N one-particle or- 
bitals whose kinetic energies are the lowest. This picture 
still qualitatively conforms to our intuition for a realis- 
tic superconductor in which the kinetic term cannot be 
viewed as a small perturbation to the pairing term. This 
suggests that insights gained in an investigation of the 
strong pairing limit gives qualitative information about 
superconducting correlations in a real superconductor. 



V. TIME-DEPENDENT PERTURBATION 
THEORY 

This section is devoted to the perturbative treatment 
of time-dependent phenomena. We consider a system 
that, in the absence of any perturbation, is described 
by the pairing Hamiltonian in the seniority limit P = 
— G5+5. At a time t = a, perturbation of the form 

j 

V= (43) 

is switched on. 

We will calculate the linear response of the expecta- 
tion value of an additive one-body operator A to this 
perturbation. From standard time-dependent perturba- 
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tion theory follows that, to first order in the perturbation 

{A) (t) = {z\ Aj{t) \^) + 1 ^* dt' {^\ [Ajit),Vj{t')] \^) , 

(44) 

where Aj{t) — exp (iPt/h) Aexp {—iPt/h) and a similar 
definition holds for Vj. Here \i) is the state the system 
was in at time t = 0, which we will take to be the un- 
pertm'bed ground state \i) oc (S'+)^|0). Note that if 
A commutes with the unperturbed Hamiltonian P, and 
\i) is a non-degenerate eigenstate of P, then it is also 
an eigenstate of Aj(t), and hence the expectation value 
of the commutator in H44|) vanishes. Thus, any additive 
one-body operator such as the current operator that is 
odd under time-reversal and therefore commutes with P 
(see the discussion following H15|l in Section ^ has no 
linear response to an external perturbation such as V. 

To compute the matrix elements involved in (|44|) . we 
will work with the mapped system. Consider for instance 

{i\Ai{t)Vi{t')\t) ^ i-\AiD{t)ViDit')\-) , (45) 

where Ajuit) = exp (iPnt/h) Ajj exp {—iPot/h) and a 
similar expression holds for Vioit')- On the right-hand 
side of H45|l . we can insert the identity operator in 
its spectral representation in the boson-fermion Hilbert 
space, between Ajd and Vjd- The operators Am and 
ViD, just like their 'Schrodinger' analogues, only have 
non-zero matrix elements between eigenstates of the un- 
perturbed Hamiltonian that lie in the same or adjacent 



levels. Hence the sum over states in the spectral repre- 
sentation of the identity operator may be restricted to 
the levels s = and s = 1. As it can easily be seen 
that the s = contributions from the two terms in the 
commutator cancel we focus on the s = 1 terms such as 

^ i-\AiDit)\i^,^.)Ui,^.\ViD{t')H. (46) 

The states |/ii^*) as defined in Section Hvl have two ideal 
fermions, one with spin pointing up and one with spin 
pointing down. This is the only possibility we need to 
include, because the operators A and V are assumed to 
be unable to cause spin flipping. This guarantees that 
their images, when acting on |— ) only create ideal ferniion 
excitations with as many spins pointing up as pointing 
down. 

Since |— ) and |/iJ^*) are eigenstates of the image of the 
unperturbed Hamiltonian Pq, we have 

(-|A,z5(i)|M^0 = exp^HAD-{t)\f^'^*),m 
i^il^,\AIDit)\-) = exp'-^{fii^,\AD+it)\-), (48) 

and similarly for the matrix elements of ViD{t'), G being 
the energy difference between the s = and s = 1 levels 
©. Thus it holds that 



{z\[Aiit),Viit')]\z) = ^2 exp'-^^^^^{-\An^\^^^^.){^l^^.\VD+\-) 

_ exp (_| vn- Ifi,^.) (a*!.. I An+ \-) . (49) 

At this point the calculation becomes a matter of substituting the correct expressions for the matrix elements and 
integrating over t' . Finally we find, to first order in V 

{A) [t) ^ + ^^^^ _ 



-,v) > si- 



sin^— . (50) 



Note again the N{Q, — N) factor in front of the first- 
order correction. As was explained in the section on time- 
independent perturbation theory, this factor accounts for 
exclusion-principle related deviations from bosonic be- 
haviour by the Cooper pairs. In the context of the per- 
fectly regular metallic nanograin this has the following 
implication: Let us suppose that the degeneracy 251 of 
the partially filled Fermi-level is large compared to the 



number N of Cooper pairs. If we increase the num- 
ber of Cooper pairs by adding electrons to the system 
while keeping the temperature fixed, the linear response 
to (say) an external electrical field initially grows, but 
not linearly. The increase in the linear response with 
an increase of particles tapers off as more electrons are 
added. When the Fermi-level is half-filled the addition of 
another Cooper pair does not change the linear response 
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of the system at all. Then, as more Cooper pairs are 
added, the linear response actually decreases. 

To see if the result (|50|) makes sense, we consider a one- 
dimensional system of 2N particles on a ring of length L 
(with periodic boundary conditions). In this case the 
index v S + l,...,j is an integer, labeling a 

one-particle wave function (j)f^{x) = -^exp{2ni^x/L). 
We will take V to be the external potential V{x) = 
Vbcos(27ra;/L), so that = Vb ((5^,^+i -f (5^+1,,.) /2. 
We will calculate the linear response of the density-of- 
particles operator 

3 

P{x)^ 0*(x)(/),.(x) (a++a,,+ -Ha+_a,._) , (51) 



so that Pfiv{x) = exp 2llillzJiEL^ One may ask if it is 
sensible to calculate the linear response of the density of 
particles. After all, the current has no linear response. 
Would the continuity equation —dp/dt = dj/dx not im- 
ply that the density should then also have a zero linear 
response to V7 It would, but the continuity equation 
simply does not hold for the pairing Hamiltonian. (Re- 
call that the continuity equation is usually derived for 
a Hamiltonian with a momentum-squared one-body ki- 
netic term and a two-body inter-particle potential that 
depends on the distance between pairs of particles only. 
The pairing Hamiltonian P looks nothing like this.) 

We consider the cut-off j to have some physical ori- 
gin, as opposed to it being the result of a regularization 
procedure. For instance, for a large or irregular grain, 
27r^j^?i^/mL^ would be set equal to the Debye energy of 
the lattice. For a small regular grain it is equal to the 
degeneracy of the partially filled Fermi-level. The cut-off 
is finite and will not be sent to infinity at the end of the 
calculation. Strictly speaking, the operators p and V do 
not leave invariant the fermion space of states with one- 
particle indices inside the cut-off only. For instance, V 
has, among others, a non-zero one-body matrix element 
Vjj+i. It is entirely possible to extend the formalism we 
have developed to accommodate this situation. However, 
this only leads to small corrections if the cut-off is large. 
We therefore simply set the offending matrix elements 
V±j,±j±i and V±j±i^±j equal to zero. 

We first compute 



{p^iv + P-u~t,){Vv^, + V-^,-^) 



2Vn ^ 2n{v-u)x 

cos 

ij 



2Vo 



L 

4Fo 



(S2 — 1) cos . 

L L 



(52) 



The third line was obtained from the second by noting 
that there are 2j = 17 — 1 terms in the summation over 
p and V, namely {p,iy) = {-j + I, -])■,{- j + 2,-j + 
1), . . . , (j, such that <5^,iy+i is not zero, and similarly 



for ^^+1 J,. Noting that tri V — 0, we therefore find that 
the general formula of (|50|l becomes 



2N NVo 



N\ 2nx oGt , . 

2^' 



in this specific instance. The spatial factor — cos {2ttx/ L) 
indicates that, as expected, particles disperse from the 
potential hill and accumulate in the potential valley. 

As another application of time-dependent perturbation 
theory, we consider a system that, before time t — 0, 
is described by the Richardson Hamiltonian, but with 
the pairing constant G so large that the one-body ki- 
netic term could be ignored. The system is assumed to 
be in the seniority limit ground state \GS) oc (5*+)^ |0) 
at i = 0. Then, at time t — 0, the pairing interac- 
tion strength is reduced to a value such that the kinetic 
term has to be taken into account. Such a tunable pair- 
ing interaction might for instance be encountered if our 
system is an ultra-cold trapped atomic Fermi gas. We 
mention in passing that recently the pairing Hamiltonian 
P = —GS~^S was discussed in precisely this context 
With the aid of the Dyson mapping the following result 
is easily obtained. The probability for the system to be 
found in the seniority limit ground state at a time t after 
the pairing strength was lowered, is to first non-trivial 
order in the kinetic energy 



Pgs^gs = 1 - 16iV. 



n - 1 



■ 2 Gt 
G^"" 2K- 



(54) 



In this expression, the by now familiar N{n — N) factor 
instead of simply TV again confirms that the likelihood 
that a Cooper pair will break up into unpaired fcrmions 
depends on the number of availiable unoccupied states for 
the electrons to scatter into, and hence on the number of 
Cooper pairs in the system. 

A matter that has not been addressed yet, is how to 
vary the number of Cooper pairs in the ground state 
without varying the number of availible states fl. While 
changing the occupancy of the degenerate Fermi level in 
a very small, regular nano-grain is feasible, the ratio N/D, 
in most other systems is not under experimental control. 
Rather, the one-particle levels involved in Cooper-pairing 
lie in a symmetric window above and below the Fermi- 
energy. The size of the window is determined by the 
electron-phonon interaction, and of the order of the De- 
bye energy, i.e. ^ lO^^cV. However, suppose one could 
find a system where the Fcrmi-cnergy is very close to the 
bottom or top of the conduction band. We will discuss 
the case of a Fermi-energy close to the bottom of the 
band. This could be attained in an electron doped semi- 
conductor. By varying the doping concentration one can 
shift the Fermi-energy up or down. Thus a situation can 
be realized where the energy-window for levels to par- 
ticipate in Cooper-pairing extends below the conduction 
band. There will then be fewer occupied levels than unoc- 
cupied levels in the energy-window, simply because there 
are no states in the section of the window that extends 
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below the conduction band. The ratio N/il will be less 
than a half, and furthermore, it is tunable: by varying 
the doping concentration, the number of electrons in the 
conduction band, and hence the Fermi-energy and the 
number of states within the relevant energy window can 
be increased or decreased. Of course, one needs to find 
a semi-conductor in which the lack of screening of the 
electron-electron interaction does not override the pair- 
ing interaction. 



one-particle levels are the same. This picture agrees 
qualitatively with our intutition also for a real su- 
perconductor in which the kinetic term cannot be 
treated as a perturbation. Since there is a smooth 
transition from the stong pairing limit to more re- 
alistic values for the pairing interaction strength, 
we expect that results such as the non-extensivity 
of the linear response also holds for the nanograins 
studied experimentally. 



VI. SUMMARY 

The salient points in this paper are the following. We 
used the Dyson mapping as a similarity transform to map 
the Richardson Hamiltonian onto an equivalent boson- 
fcrmion Hamiltonian. This Hamiltonian uses bosons to 
describe the bound fermion pairs that result from the 
pairing interaction, instead of the time-reversed collec- 
tive fermion pairs created by the operator . From the 
structure of the mapped operator emerges our first main 
result: An additive one-body operator only has non-zero 
matrix elements between eigenstates of the seniority pair- 
ing model, if these eigenstates come from the same or 
adjacent energy levels in the spectrum of the pairing 
Hamiltonian. With this selection rule it becomes feasible 
to investigate the strong pairing limit of the Richard- 
son Hamiltonian without relying on the Bethe-Ansatz 
solvability of the Richardson Hamiltonian. While our 
aim remains the same as that of previous authors, i.e. 
to develop a simple and intuitive analytical picture of 
the low-temperature physics of superconducting metallic 
nano-grains, the new tools at our disposal allow us to 
obtain new results: 

1. In Section Hvl we calculated the ground state ex- 
pectation value of an arbitrary one-body operator 
for a system of 2N electrons in the strong pairing 
limit. In Section El we showed that the selection 
rule we identified and used in conjunction with the 
Dyson mapping also facilitates time-dependent per- 
turbation theory. Here we explicitly perturbed the 
strong pairing limit with a more general perturba- 
tion than the kinetic term in the Richardson Hamil- 
tonian. Results obtained in both calculations show 
that the likelihood for a perturbation to break up 
a Cooper pair depends on the number density of 
Cooper pairs in the system. This is explained by 
the fact that the electrons that constitute a Cooper 
pair obey the exclusion principle and underlines 
that there are circumstances in which Cooper pairs 
do not behave as bosons. 

2. Also in Section ^3 we calculated the occupation 
probability of single particle levels in the strong 
pairing limit. We found that a picture emerges of a 
competition between the kinetic term that favours 
a Fermi-sea, and the pairing term that favours a 
situation where the occupation probabilities of all 
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APPENDIX: DEGENERATE PERTURBATION 
THEORY 

In this appendix the general formalism is presented in 
detail, for finding perturbatively the spectra of (possibly 
non-hermitian) operators if there is degeneracy when the 
perturbation is switched off. The argument is developed 
along the lines of the non-degenerate (hermitian) case as 
treated by Sakurai,^'' We could not locate any source that 
gives a satisfactory treatment of the degenerate case. 

We are given an operator H{X) ~ Hq + XHi. Neither 
Hq nor Hi have to be hermitian. We assume that H is 
fully diagonalizable, at least in a finite region around A = 
0. For A ^ 0, we assume that H{X) has a non-degenerate 
spectrum, so that a single label a = 1, 2, . . . , M may be 
used to specify eigenvalues and eigenvectors uniquely: 

H{X}\an)^ = E^{X)\an)^, (A.l) 

with Ea{X) + Ep{X) if a 7^ /3 and A 7^ 0. The sub- 
script R is used to indicate that we are dealing with right- 
eigenvectors, later to be contrasted with left-eigenvectors, 
that do not coincide with the right eigenvectors when 
H{X) is not hermitian. We take it that all the Ea{X) 
and are analytical functions of A in a finite region 

around A = 0. Furthermore, if two eigenvalues converge 
at A = 0, i.e. Ea{0) — £^/3(0), we assume that these 
eigenvalues have different first derivatives at A = 0. This 
condition is referred to by saying that all degeneracy is 
lifted in the first order. At A = 0, the spectrum of H may 
be degenerate. Our notation takes care of this as follows: 
At A = we let M' be the number of distinct eigenval- 
ues of Hq. We then indicate these distinct eigenvalues 
of Hq as £i, I = 1,2,...,M'. We partition the labels 
a — 1,2, . . . , M into disjoint sets Pi, labeled by integers 
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1 = 1,2,... M' , by defining Pi to be the set of all indices 
a such that Ea flows to £i as A goes to zero. We can then 
also define a function p from the index set 1, 2, . . . , M to 
the index set 1,2,..., Af as follows: p{a) = / where I 
is the unique label such that a £ Pi. We refer to the 
space spanned by all eigenvectors of Hq with eigenvalue 
£i as the £i eigenspace of Hq. It is assumed we know a 
basis for each £i eigenspace of Hq. The basis elements 
that span the £i eigenspace will be denoted 1-0, ^_r) where 
^p is a discrete index that runs from 1 to some integer 
oJi which is the dimension of the £i eigenspace. The set 
{Itp, In) : I = 1, . . . , M'; = 1, . . ■ , wz} forms a basis for 
the domain of H{X). It is in terms of these basis states 
that we want to express the eigenstates of H{\) approxi- 
mately. Corresponding to this basis we can always find a 
left basis {I?/", Il) : / = 1, . . . , A/'; -0 = 1,..., uoi} such that 
the equations 



(A.2) 



hold for all I, m, ip and 4>. 

We uniquely define a projection operator H; as the lin- 
ear operator that maps any vector in the £i eigenspace 
of Hq onto itself, while mapping any vector that can 
be written as a linear combination of vectors from the 
other eigenspaces of Hq, onto zero. (This might not be 
an orthogonal projection operator though. When Hq is 
non-hermitian, there are £i eigenspaces that are not or- 
thogonal to each other.) An operator 11; is defined as 
Hi = I — 11;. In terms of the states and |0, mj^), 

the identity operator I, the unperturbed operator Hq and 
the projection operators 11; and 11; can be expressed as 
follows: 

M' uji 

1=1 ^=1 

Ho - ^£;^|V,^i?)(V',ZL|, (A.4) 

n, = ^|v^,z,i)(v^,ZL|, (A.5) 

M' uik 
Note that IljIIm = SlmJ^m- 

With the above preliminaries out of the way we now 
expand the eigenvalues Ea{X) and eigenvectors |Q;ij)^ of 
H{\) in terms of A: 



Ec.{\) = + AAW + A^Ai^) + . . . 



,.(1) 



(2) 



(A.7) 



This definition fixes the directions of the A-independent 



zero'th order eigenvectors 



a 



(0) 



uniquely as the direc- 



A goes to zero. This direction is well-defined since we 
assumed that H{X) has a non-degenerate spectrum for 
A 7^ 0. The magnitudes of the eigenvectors are left arbi- 
trary. 

The set of eigenvectors | a^^''^ | also forms a basis 
for the domain of H{X). We use it to define a second basis 

(sometimes called the contra- variant basis) | a^''') | 
through the set of equations 



(A., 



Hq = £ 



p(a) {a 



als 



Equation l|A.8(l implies that 

Using the states | Q!j^''^| and function, 

{ ^^^^^ } define another set of projection oper- 

ators that will come in handy during our calculations. 
Define a projection operator 



JO) 



JO) 



(A.9) 



that singles out the 



a 



(0) 



component of the decomposi- 

M 



tion of any vector in the | ct^"*^ | basis. Let tTq, be 

the complementary projection operator of tTq in the £p(a) 
subspace by setting 



/3ePp(„,\{a} 



(^: 



= n 



p(a) 



(A.IO) 



Thus, the identity operator may be decomposed as 



I = n 



p(a) 



(A.ll) 



The next step is to choose the relative normalization be- 
tween the eigenstate |afl);^^ and Icti);;,,. The most conve- 
nient choice turns out to be {^af^^ ajij =1. With this 
choice it follows that 



JO) 



(A.12) 



Now we rewrite the eigenvalue equation in the form 

- Hq) \an),^ = (Ai/i - A„(A)) \aR)^ , (A.13) 

where A„(A) = Er=i ^i'^A^ 

The first thing we do with this equation is multiply it 
by np(Q), noting that np(„) {£p{a) - Hq) = 0, to find 



np(„) (Affi-A„(A))|aflj^ = 0, 



(A. 14) 



which has to hold order for order in A. Note that, by 
definition, the zero'th order term in the A-expansion 



tion to which the a-eigenvector of -ff(A) converges as of |aij)^, namely 



„(o) 



is an element of the £, 



p{a) 
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.(0) 



JO) 



Keep- 



eigenspace of Hq and hence llp(a) 

ing this in mind, we look at the terms in (|A.14|I of lowest 
order in A to find 



a 



(0) 



All) 



„(o) 



(A.15) 



This eigenvalue equation enables us to find the zero'th 
order vectors a^"*^ as the eigenvectors of Ilp(^a)Hi in the 
^p{a) eigenspace. Since we assume all degeneracy to be 
lifted in the first order, no two eigenvalues Aa^ and A^^-* 
with a, (3 £ Pp(a) are the same, and (IA.15P is necessary 



and sufficient to determine the various 



a 



in the £, 



eigenspace. 

In order to proceed, we define an operator 



1 



— ^ Oh 

l^k " 



p{a) 



(A.16) 



We return to ljA.13p and multiply it with Qp(c), noting 

that Qp(a) {£p{a) " Hq) = llp(a)- Thus wc arrive at one 
of the most important formulas in the present discussion, 
namely 

np(a) \aR)x = Qpic) (Ai?i - A„(A)) \aR)^ . (A.17) 

This equation has, on the left-hand side, the part of the 
eigenvector that can be written as a linear com- 

bination of vectors in Si eigenspaces with I ^ p{a). We 
want a similar equation for that part of that lies 

inside the £'p(a) eigenspace. 

For this purpose, we firstly rewrite l|A.14ll to read 

(Aq(A) - \Iip(^a)HlI\-p{a)) ^p{a) \(^r)\ 



Anp(„)i7inp(„) \aR)^ = 0, (A.18) 



by recalling that Tlp {a) + ^ p{a) — I and {jlp(a)) = ^p{a 



Now we substitute (IA.17I1 into the second term of lA.l 
to find 



(Aq(A) - \Iipi^a)Hiilp(a)) np(Q) \a.R)^ 

- Xnp(^c.)HiQp(c.) {XHi - A„(A)) \aR)^ = 0. 



(A.19) 



Define an operator 



/3e-Pp(c)\{"} 



^(1) 



.(1) 



(A.20) 



Note that qa ^Aq ^ — Tlp(a)HiIlp(^a)j = ^^a and that 
Qa^p{a) = Qa- Hcncc, by multiplying ljA.19|) with Qa 
and dividing by A, note that 



qaHiQp^a) (XHi - Aq(A)) \aR)^ 

OO 

-Y,Ai'+'h\^\aR),. (A.21) 



fc=i 



This is the counterpart of (|A.17|I that we needed for the 
part of lai^)^^^ that lies inside the £ p(a) ei genspace. When 
we combine the results of (|A.12|I . IjA.lTp and (|A.21|I . re- 
membering that the identity operator may be decom- 
posed as I = Ilp(^a) + tTq; + tTq, we find 



I««)a = "fl^j + (1 + Qpia) (Affi - A„(A)) la^);, - ^ A^'^+^U'^g^ \aR), 

Collecting terms of order iV -|- 1 we find, for > 0, 



fe=i 



a 



(JV+l) 
R 



(1 + qaHi) Qpla) { Hi 



p[c 
N+2 

M=2 



AN) 



(N+2-M) 



r 



Af-l-l 
M=l 



{N+ 



1-M)^ 



(A.22) 



(A.23) 



Note that Qp(a) (^^r' j = and qa 
for A'^ > 1, the upper bounds on the two summations 
in the above expression may respectively be decreased 
from A -I- 1 to A^ and from A^ -I- 2 to A^ -I- 1, while for 
the case where A^ = 0, the two summations may be left 



„(") 



0. Thus, 



out entirely. To complete the expansions, we need an 
expression for a[^^ . This we find by multiplying IjA.lSp 

, recalling our normalization 



with the functional ( a 



convention ( a^^^ 



4°^ 



aR ) =1, and collecting terms of order 



16 



M in A to find 



Our final results are then: 



(0) 



(A.24) 



a 



(w+i) 

R 



(1 + qaBx) Qpi^a) 
N+l 



M=2 



(N+2-M) 



N 
M=l 



^R 



(A.25) 



r 



and 



Hi 



a 



(N) 
R 



(A.26) 

and 
exceed 

the upper bounds, the summation contains no terms, 
i.e. Qfo'' 1 = + QcyH^) Q„(n\Hi ao"* )• The formulas 



It is understood that, in the case where N 
the lower bounds of the summations in (IA.2 



(0) 

of ((223 and (|X26)l are all we need to calculate the 
corrections of order + 1 to both the eigenvectors and 
eigenvalues of H, if the corrections of orders 0, 1, 2, . . . , A^ 
are known. Together with the initial conditions provided 
by the eigenvalue equation 



a 



(0) 



v(o) 



(A.27) 



this allows us to write down in principal the expansion of 
any eigenvector and its eigenvalue to arbitrary order in A. 
Concerning ljA.26(l for the eigenvalue-correction of order 



A^+l, note that ( a 



JO) 



HiQa = 0. This is true because the 



states 



JO) 



and functionals 



with a, (3 Pi were 



sense that 



Hi 



chosen to diagonalize Hi in the £'p(c,) eigenspace, in the 
,.(0)^ —Aa^Sa.p- This means that, 

("^^^ in ^rm . the 



when we substitute from (|A.25|1 for 

terms that are multiplied from the left with an operator 
Qa disappear and we are left with 



(0) 



AiN+2) ^ 



N 



M=l 



HlQp[a)Hi 

HlQp(a) 



a 



(N) 



„(0) 



a 



(Af+l-M) 
R 



(A.28) 



which is an expression for the eigenvalue-correction, with 
all the dead-wood cut away. If an unperturbed eigenstate 

q;^'^ has a non-degenerate eigenvalue, i.e. Hq does not 

contain other eigenvectors with the same eigenvalue as 

ifp'' ) , the situation and hence the formulas are simpler 



than above. In this case, the proper operator to use for 
Qa is simply the zero-operator. With this adjustment the 
general formulas reduce to formulas valid for the specific 
case of no degeneracy. 

Below, we write down explicitly the first order correc- 
tion to an eigenstate and the second order correction to 
its eigenvalue, for the general case where degeneracy is 
present in the unperturbed system: 



E 



UiHi 



--p(a) 



■ E E 





HiIliHi 













E 

l^p(a) 



„(0) 



HiUiHi 



,,(0) 



p{a) 



£i 



(A.29) 



(A.30) 



Note that the formula for the second order correction to 



vector 

n. 



the eigenvalue, (|A.30p looks the same for the degenerate 
case as for the non-degenerate case. The only difference 
is that, in the degenerate case, the zero'th order eigen- 

a^"*^ is a simultaneous eigenvector of Hq and 

p{a)Hi, lying in the £p(^a) eigenspace, whereas in the 
non-degenerate case it is uniquely fixed by diagonalizing 
Ho- 

When we apply this perturbation expansion method 
in conjunction with the Dyson mapping, there is still the 
matter of the physical subspace to address. The question 
is: If we approximate the eigenvectors of H{\) by includ- 
ing only a few terms in the perturbation expansion, do 
we lose the notion of physical eigenstates that lie exactly 
in the physical subspace? Recall that if we diagonalize 
H{X) exactly, we can express the physical subspace as 
the space spanned by a certain subset of the eigenvec- 
tors of H{\). If we only approximate the eigenvectors 
by the first few terms of their A expansions, can we still 
find a subset of these approximate eigenvectors that span 
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the physical subspace? Or do the A-expansions for the 
physical eigenvectors, if truncated after a few terms, only 
lie close to, but not necessarily in the physical subspace? 
Perhaps not unexpectedly, the answer is that the approx- 
imate physical eigenvectors still lie exactly in the physi- 
cal subspace. This can be seen as follows. Any operator 
that leaves the physical subspace invariant is said to be 
a physical operator. Then if we assume that H{X) is the 
linear extension of the image of a fermion operator for all 
A, it follows that Hq and Hi are both physical operators. 
Furthermore, it is a simple matter to show that the oper- 
ators 11; and Qi are physical operators. This means that 
HiHi is also a physical operator. We find the zero'th 



order eigenvectors 



by diagonalizing the operator 



Ilp[a)Hi in the £p(^a) eigenspace of Hq. Since IIp(^a)Hi is 
a physical operator with a non-degenerate spectrum, the 
eigenvectors thus obtained can be separated into a set 
that spans the overlap of physical subspace and the Sp^a) 



eigenspace, and a set of ghost states. Since the operator 
Qa is defined in such a way that it has the same eigen- 
vectors in the £p(^a) eigenspace as the physical operator 
np(a)-ffi, it follows that Qa is a physical operator. Then, 
if we start the recursion of (|A.16|I with a state from the 
physical sector, the corrections that the recursion gen- 
erates remain inside the physical subspace. This follows 
inductively from the fact that higher order corrections 
in ljA.25|) are generated by acting with physical opera- 
tors on lower corrections. If the zero'th order eigenvector 
is physical, this implies that the first order correction is 
physical. If the zero'th order eigenvector and the first or- 
der correction are both physical, then so too is the second 
order correction, etc. This implies that the exact phys- 
ical eigenstates, when expanded in terms of A, consist 
of terms that all lie in the physical subspace themselves. 
Truncating the expansions after a finite number of terms 
therefore still leaves one in the physical subspace. 
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This becomes an issue when the finite set of one-particle 
indices ^ = —j, ■ . ■ ,j is a subset of a larger (possibly in- 
finite) set of indices. See for instance the example of the 
particles confined to a ring, discussed in Section IVI 



